Effective theory of rotationally faulted multilayer graphene - the local limit 
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Interlayer coupUng in rotationally faulted graphene multilayers breaks the local sublattice- 
symmetry of the individual layers. Earlier we have presented a theory of this mechanism, which 
reduces to an effective Dirac model with space-dependent mass in an important limit. It thus makes 
a wealth of existing knowledge available for the study of rotationally faulted graphene multilayers. 
Agreement of this theory with a recent experiment in a strong magnetic field was demonstrated. 
Here we explore some of the predictions of this theory for the system in zero magnetic field at large 
interlayer bias, when it becomes local in space. We use that theory to illuminate the physics of 
localization and velocity renormalization in twisted graphene bilayers. 
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I. INTRODUCTION 

Experiments indicate that the 10-100 individual 
graphene layers grown on the carbon-terminated 
face of SiC are surprisingly well decoupled from 
one another electronically. Early spectroscopic 
measurement^!^ found a linear low-energy electronic 
dispersion to the experimental precision, like that of 
single-layer graphene^i^. In scanning tunneling mi- 
croscopy /spectroscopy (STM/STS) measurements the 
Landau level quantization of the material in a mag- 
netic field was found to be essentially that of single-layer 
graphene^. Theoretically it has been shown that this 
approximate decoupling of different layers is due to a rel- 
ative twist of the layers with respect to each other— —. 
A renormalization of the electron velocit}^!^^, van Hove 
singularities, and interlayer transport^i^ have been dis- 
cussed as residual effects of the interlayer coupling. 

In a recent STM measurement on multilayer epitaxial 
graphene^^ a spatially modulated splitting A ^ 10 meV 
of the zeroth Landau level (LLq) was observed. In view 
of the above this finding is intriguing, since the states 
forming LLq of an isolated layer of graphene without 
electron-electron interactions are degenerate. Many as- 
pects of the experimental data indicate that this split- 
ting is due to the coupling between graphene layers. In 
Ref. [m we proposed a phenomenological theory of the 
interlayer interaction. In that theory a "staggered" elec- 
tric potential (a potential with opposite sign on the two 
sublattices) breaks the sublattice-symmetry locally. This 
model qualitatively accounts for the main features of the 
experimental data. In Ref. we have presented a mi- 
croscopic theory of the interlayer coupling in rotationally 
faulted graphene multilayers that reduces to the phe- 
nomenological model of Ref. 15 in certain limits. The 
theory is formulated for a single layer of graphene and 
it accounts for the coupling to other layers by effective 
potentials and an effective mass that are possibly non- 
local in space. The theory of Ref. [H accounts for the 
main features of the experimental finding^, both qual- 
itatively and quantitatively. 

A number of intriguing results have been obtained 



theoretically in electronic structure calculations of 
rotationally-faulted multilayer graphene also in zero mag- 
netic fieldiii^. One may therefore ask whether the the- 
ory of Ref. [H can provide an intuitive understanding also 
of these results, as it did for the physics of the material in 
high magnetic field: is that theory an advantageous start- 
ing point to exploring the physics of rotationally-faulted 
multilayer graphene also in zero magnetic field? 

In this article we give a partial answer to th at q uestion 
by exploring predictions of the theory of Ref. [l6| in zero 
magnetic field for quantities that have displayed interest- 
ing features in the calculations of Ref. |Tl| and by seeking 
an interpretation of the results in qualitative terms. We 
focus on the spatially local limit of the theory uM that 
corresponds to the phenomenological model of Ref. |l5|: a 
single-layer Dirac model with oscillating effective poten- 
tials and a space-dependent mass. That limit is realized 
in the presence of a large interlayer bias. The theory pre- 
dicts a density of states in qualitative agreement with ex- 
perimental topographic STM measurements. Moreover, 
our calculation qualitatively reproduces some of the main 
observations of the mentioned electronic structure calcu- 
lations of twisted graphene bilayers^ ^ such as a local- 
ization of electronic states and a corresponding velocity 
suppression. The agreement is not quantitative, since the 
calculations of Ref. ^ were not obtained in the spatially 
local limit assumed here. But in the framework of the 
theory of Ref. [H these predictions do have an intuitive 
explanation in terms of known results about the Dirac 
equation with a space-dependent mass. This suggests 
that this theory is indeed an advantageous starting point 
for the exploration of the physics of rotationally faulted 
graphene multilayers also in zero magnetic field. 

We start our discussion with Section |TT1 where we re- 
state the model on which our earlier theory^^ is based. 
In Section irm we take the limit of a large interlayer bias, 
when the effective theory of Ref. [ifi becomes local in 
space. We then proceed to evaluate the density of states 
and the electron velocity renormalization predicted by 
this theory in zero magnetic field. We do that first per- 
turbatively in the interlayer coupling in Section |lVl In 
Section |V] we then analyze nonperturbatively a toy model 



that resembles the original theory, but assumes a simpli- 
fied spatial structure of the effective staggered potential. 
We conclude in Section IVTl 

II. MODEL 

In this Section we recall the model of Ref. 16, which 
underlies also the present article. We analyze the elec- 
tron dynamics in a graphene layer "0" when coupled 
to a second layer "1," twisted by a relative angle 
{0 = 0° for aligned honeycomb lattices, cf. Fig. [T]), ne- 
glecting electron-electron interactions. The correspond- 
ing dynamics in multilayers at perturbatively weak in- 
ter layer coupling, such as in the experimenlji^, are ob- 
tained by summation over all layers coupled to the top 
layer 0. Twisted graphene bilayers have been described 
beforeiiSri^ by a tight-binding model with a local inter- 
layer coupling Hamiltonian that has parameters fitted to 
experimentpi^, 

i^int = j dr^^^^\r)T{r)^f^^\r)^h.c. (1) 

Here, the spinors ^^^^ collect the amplitudes for electrons 
on the two sublattices of layer j G {0, 1}. The inter layer 
coupling r has contributions at wavevectors h^^^ — h^^\ 
where h^-^^ are reciprocal vectors of the graphene lat- 
tice in layer j^. The Fourier components of F quickly 
decay with increasing wavevect o r ^i ^ Q i ^ ^ . We therefore 
neglect all but the zero wavevector component, setting 
F(r) = 7. In the "first star approximation" of the 
wavefunctions employed below, the distinction between 
commensurate and incommensurate interlayer rotations 
then disappears. This approximation is valid for ener- 
gies £ ^ where V is set by the Fourier components 
of F that directly connect K-points of the two layer^. 
We take the hmit < 6> <C 1, when V <C 7 (in the 
experimentfi^ ^ 0.25° and according to the estimate 
V ^ O'^'y of Ref. [H this approximation is justified at all 
accessed energies). 

In our limit < ^ <C 1 a long-wavelength description 
is appropriate, where the isolated layers j are described 
by Dirac model Hamiltonians (we set h = 1) 

H^^^ =v f dr^V^^^')^(r) [a^ ■ {-iV ^ eA{r))]iljl^\r). 

(2) 

Here, cr^y = {vcFx^ a^y) is a vector of Pauli matrices, = ± 
is the valley spin, — e the electron charge, and v the elec- 
tron velocity in graphene. We have included an external 
vector potential A to describe a perpendicular magnetic 
field B. Eq. ([2j) acts on the long- wavelength spinors 
defined by ^^f^(r) = E^. ^M^^(^)V^Mi^(^)- We write the 
Bloch functions u^l{r) = {Ep exp[zK^i) • (r-r[f ^)]}/V3 
in the "first star approximation" appropriate for the in- 
terlayer coupling problem^^. Here, p sums over the three 
equivalent Brillouin zone corners K[^^ that form the 
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FIG. 1: (color online) Moire pattern created by two graphene 
lattices with a relative twist. Top layer A/B sublattice 
atoms are shown as small blue/cyan (dark/light) spheres 
and connectors; bottom layer A/B atoms are shown as large 
red/yellow (dark/light) spheres. A region of A A alignment 
lies at the center, where each top-layer atom has a neighbor 
in the bottom layer. The A A region is surrounded by three 
AB- and three BA-aligned regions where atoms on only one 
top-layer sublattice have direct neighbors in the bottom layer. 
As a consequence, the sublattice-symmetry is broken locally. 

Dirac point of valley and r^^^ gives the position of 
an atom on sublattice ji G {A, B} within the unit cell of 
layer j. In the long- wavelength theory (which neglects 
inter- valley processes) the interlayer coupling reads 

^fint = / dr^^/>(o)t(r)t.(r)V(^)(r) + h.c, (3) 

with a matrix t whose long- wavelength components have 
wavevectors bK^^ = {Rq — 1)K^J . Here, Re is a rotation 
around the z-axis by angle 0. Retaining only those long- 
wavelength parts of t we find 

t/;A.'(^) = 2 ^e'^^--^+^^--K°'-^t' ), (4) 

o 

p 

where terms of order are neglected, while terms of order 
6Kr are kept as they may grow large. 

III. EFFECTIVE THEORY 

We next integrate out layer j = 1 in order to arrive at 
an effective Hamiltonian Hq^{uj) = Hq -\-SHq^{uj) for the 
top layer j = 0, with 

5Hf{^) = H,^t{oo + F - i/i)-^i/int. (5) 

We include an interlayer bias V that accounts for differ- 
ent doping levels of the two layer In general, Hq^ is 
nonlocal in space and it depends on the energy oo. In the 
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limit of a large interlayer bias, however, \ V\ u^j^Ov/a^ 
the sum uj -\- V — Hi becomes momentum- and energy- 
independent to a good approximation. The spatial non- 
locality and the energy-dependence of H^^ then may be 
neglected and Hq^ becomes a conventional Dirac Hamil- 
tonian (j2j) with a matrix potential 

SHf = [ d.^^(o)t(,)M£l!^^(o)(,), (6) 

which we parametrize as 

U{r)tt{r) ^ v^^{r) + i.vea,-A'^''{v) + m''^{rya.. (7) 

The interlayer coupling in this limit generates effective 
scalar and vector potentials V^^ and A®^, respectively, 
and a mass term ex a^m^^'U^ that implies an effec- 
tive staggered potential = rrf^v^ in locally Bernal 
stacked regions. It follows from Eq. (j4|) that SH^^ oscil- 
lates in space with wavevectors G = {Rq — where 
b is in the "first star" of reciprocal lattice vectors of 
graphene. We plot 5Hq^ in the parameterization of Eq. 



([7j) in Fig. O The effective Hamiltonian (|6]) promises rich 
physics. In particular the effective mass term is expected 
to have profound imphcations such as topologically con- 
fined states^^'^^. In Ref. 16 we have shown that the above 
effective theory qualitatively and quantitatively accounts 
for many features of the experiments^, which was done 
in a strong magnetic field. In the remainder of this arti- 
cle we explore some of the consequences of the effective 
potentials ([7j) in zero magnetic field. 



IV. PERTURBATIVE RESULTS 

We first explore the perturbative limit of weak inter- 
layer coupling 7^ <C Vv5K with correspondingly weak 
effective potentials, Eq. ([7]). To this end we do pertur- 
bation theory in SHq^ . We first obtain the perturbative 
corrections to the low-energy density of states po(r) = 
lim^^o as probed in STM measurements. The 

lowest order correction to po(r) = lim^^Q p{t, e)/e van- 
ishes. The leading spatially varying contribution appears 
at second order in SH^^: 



^/>o(r) 



1 



d(pk 



(5, k, v\mf,\s\ G^ v){s\ G^ ^|r)(r|5^ G, v){s' , G, v\mf,\s, k, v) 



+2Re- 



^2|G||G'| 

5, k, v\y){y\s\ G - G^ G - G^ uISHS^W, G, G, iy\SHS^\s, k, jy) 

v^\G\\G-G'\ 



Here, |s, k, i/) is an eigenstate of Eq. (|2]) at A = with wavevector k and energy sv\]i\ in valley 

1 ik r / 1 



(r|5,k,i/) = 



V2' 



(8) 



(9) 



where cpk = arctan(/c^//caj). The sums over wavevectors G, G' in Eq. ([8j) runs over all wavevectors contributing to 
SHq^. Eq. dHI evaluated for the effective Hamiltonian (|6]) results in 

eff.r^/M 1 + (|G| + |G'|)/|G-G'| 



27ru4|G||G' 



We plot the resulting relative correction to the density of has21 
states 6po/po = 27rv'^dpo in Fig. [3l The result compares 

well with the typical moire patterns observed in STM _ d k, i^|(5i^oS^|s', k — G, 

topography. This suggests that density of states correc- ^ |k|^o d\k\ ^ 'i;(|k| — s'|k — G|) ' 

tions due to the effective potentials Eq. ([7j) may be one ^ ' 

of the mechanisms that generate these patterns, besides 
simple geometric height variations of the top graphene 
layer (which would be the most straight forward inter- 
pretation of topographic STM maps). . _ / |F^^(G)p[|Gp - (k • G) 

^^k - ^2|G|4 

We next evaluate the perturbative correction to the 
electron velocity at the Dirac point in direction of the |m®^(G)p(k • G) 

momentum 'Uj^ = lim/c^o v(k) • k, where k = k/ 1 k|. One ^ 'U^IGI^ 



which, for our evaluates to 



}• (12) 
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FIG. 2: (a) Effective potential , (b) effective mass nf^, 
(c) A'f , and (d) of Eq. ^ as functions of rO/a in grey- 
scale. Scale bars span a unity increment in r^/a, where a is 
the C-C bond length (0.142 nm). Note the expected sixfold 
and threefold symmetries of V^^ and m^^, respectively. A^^ 
transforms as a vector under rotations. 




FIG. 3: Perturbative correction to the low-energy density of 
states (Spo/po) due to SHq^ . The scale bar corresponds to 
one unit in rO/a, where is the rotation angle between layers 
and a is the C-C bond length (0.142 nm). For a rotation angle 
of 3°, interlayer coupling of 7 = 300 meV, and an interlayer 
bias oi V = 400 meV (stretching the limits of our locality 
assumption) the image corresponds to a 16 nm x 16 nm area 
with Spo/po ranging from -0.125 (black) to 0.125 (white). 



Here, we have used that G • A(G) = 0, which holds in 
our approximations. We note that the effective mass sup- 
presses the electron velocity as was found for an oscillat- 
ing scalar potential V^(r) in Ref. [2l|. However, differently 
from a scalar potential, that velocity suppression is in 
the case of a mass not perpendicular to the direction G 
along which m®^ varies, but along that direction. For the 
effective Hamiltonian §^ [with b = in Eq. (j4])] Eq. ([I2j) 
evaluates to 



(13) 



which is isotropic in space. We anticipate anisotropic 
contributions to Sv^ at higher orders of perturbation 
theory. The interlayer coupling reduces the velocity, in 
agreement with earlier calculations^' for twisted bilay- 
ers at V = 0. We conclude that in the perturbative 
regime of weak interlayer coupling the predictions of our 
theory are consistent with earlier experimental and the- 
oretical work. They moreover have a straightforward in- 
terpretation in terms of previous results for electrons in 
a superlattice potentiaP^. 



V. NONPERTURBATIVE RESULTS 

We now turn to the more challenging, but also more 
interesting nonperturbative limit. For strongly cou- 
pled, twisted graphene bilayers a number of intrigu- 
ing results have been obtained in electronic structure 
calculation o^^i^^ . For instance, a localization of the elec- 
tronic wavefunctions in locally AA-stacked regions of the 
sample and a severe suppression of the electron veloc- 
ity at small twist angles have been observedii. Here, 
we show that the mentioned phenomena find an intu- 
itive interpretation in terms of the oscillating mass in a 
toy model of our effective theory. In this toy model we 
assume translational invariance in one space direction. 
Our calculation extends earlier theory of Dirac electrons 
with a (scalar) superlattice potential^^Ti^i to the case of 
a periodic mass term. 

Our toy model has a mass term, which oscillates in 
only one space-direction: 



m®^(r) = msgn [cos{x/l)]. 



(14) 



Although not directly applicable to the problem of the 
interlayer coupling in multilayer graphene, such a model 
is expected to capture some of the same physics. We 
add a term m^^v'^az to Eq. ([2|) at A = and find the 
electronic spectrum along the lines of Ref. The re- 
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FIG. 4: Bandstructure of a Dirac model subject to a mass 
term Eq. fTl)) oscillating on a length scale /. In the plot, to 
every momentum ky, the energies e where electronic states 
exist are marked white. For the plot we chose m — v/l 



suiting band structure as a function of the momentum 
ky in the direction with translational invariance is plot- 
ted in Fig. m In presence of the periodic mass the elec- 
tronic spectrum breaks up into minibands, as for periodic 
scalar potentials'^' Differently from the scalar case, 
however, we find that a periodic mass does not gener- 
ate new Dirac points, even in the nonperturbative limit 
vml ^ 1. Instead, the bandwidth of the lowest energy 
band becomes exponentially suppressed in vml with cor- 
respondingly suppressed electron velocity: the velocity in 
x-direction at zero energy is exponentially small in vml^ 



vml 



\J2 cosh'um/ — 2 



(15) 



The wavefunctions at large vml have dominant weight 
around the locations x = 27rn (with integer n) where 
m®^ changes sign. Those are the locations, where for 
an isolated kink of m®^, i.e. a point at which m®^ 
changes sign, topologically protected zero energy states 
are expected^ i'^i'^i'^ . The lowest energy band observed 
in Fig. m may be thought of emerging from hybridiza- 
tion of those zero energy states. The larger vml the less 
overlap occurs between the states localized at adjacent 
kinks of m®^. This qualitatively explains the exponen- 
tially small bandwidth of the lowest energy band in Fig. 

El 

The above findings closely resemble the above- 
mentioned observations of earlier studies of twisted 
graphene bilayers^^. Also there the wavefunctions were 
reported to be localized in the AA-stacked regions of the 
moire pattern when ^ <C 1, which implies large with 



vml 1. Those are indeed the regions, where m®^ in 
H^^, Eq. d?]), chan ges sign. They thus directly corre- 
spond to the regions where the wavefunctions in our toy 
model Eq. (p!^ are concentrated. Moreover, in Ref. |Tl| 
the electron velocity was found to be strongly suppressed 
at small ^, corresponding to large vml. Also this is in 
qualitative agreement with the prediction of an exponen- 
tial suppression of Vx by our model [Eq. (p!5]) ]. 

Although our theory is not directly applicable to the 
calculation of Ref. [Ill, because that calculation was done 
at V = and the moire pattern was naturally two- 
dimensional, it indeed appears to capture some of its 
essential physics. The above calculation thus suggests 
an intuitive interpretation of some of the prominent non- 
perturbative effects in twisted graphene bilayers in terms 
of zero energy states that are induced by the topology of 
the mass term in our model. 



VI. CONCLUSIONS 

In this article we have discussed some of the implica- 
tions of the effective theory of rotationally-faulted mul- 
tilayer graphene that was put forward in Ref. [l6|. We 
have focused on its local limit of a large interlayer bias, 
when this effective theory takes the form of a conven- 
tional Dirac model with space-dependent potentials and 
mass. While we discussed the implications of that theory 
for graphene multilayers in a magnetic field in Ref. [l6|, 
here we have explored its consequences in zero magnetic 
field. In the perturbative limit of weak interlayer cou- 
pling we found corrections to the density of states that 
are consistent with the typical moire patterns observed 
in topographic STM measurements. This suggests that 
these patterns may not only arise because of height fluc- 
tuations, but may at least partially be due to density 
of states variations. We moreover have found a velocity 
correction consistent with earlier calculations in different 
limits. 

To access the most interesting nonperturbative regime 
of strong interlayer coupling we have analyzed a toy 
model that captures most of the essential ingredients of 
our effective theory. We have demonstrated that that 
model predicts almost localized electronic states and an 
exponential velocity suppression. These predictions give 
an intuitive interpretation to prior electronic structure 
calculations for twisted graphene bilayers in terms of 
topologically protected zero energy states localized at 
kinks of an oscillatory Dirac mass term. Partially answer- 
ing the question we raised at the outset, the presented 
calculations lead us to the following conclusion: while 
our real-space theory of the interlayer coupling certainly 
is an advantageous description of twisted graphene bilay- 
ers with a large interlayer bias, it qualitatively captures 
much of their essential physics even without such bias, 
when it does not strictly apply. The theory of Ref. [l6| 
indeed appears to be an advantageous approach to the 
physics of twisted graphene bilayers, also in zero mag- 
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netic field. 
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